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INTRODUCTION 


The  goal  of  this  work  is  the  development  of  techniques  that  permit  timely  and  efficient 
evaluation  of  different  propagation  models.  The  resulting  program  is  to  be  as  sophisticated  as 
possible  with  concerns  for  computer  run  time  to  be  subordinate  to  the  use  of  realistic  models  of 
the  environment.  Employed  with  new  databases  of  high-latitude  propagation  measurements,  this 
model  will  be  used  to  evaluate  and  improve  the  faster-running  models  and  point  to  deficiencies 
in  the  ionospheric  models  and  the  high-frequency  (HF)  prediction  models.  In  addition,  the 
development  of  new  graphical  displays  will  enhance  understanding  of  the  propagation  environ¬ 
ment.  The  emphasis  of  this  development  is  on  the  personal  computer  for  portability  and  ease  of 
access.  The  current  operating  system  is  OS/2,  Version  2.1.  One  of  the  modes  of  graphical  output 
is  specific  to  this  operating  system. 

SPATIAL  SMOOTHING 

The  first  step  in  the  development  of  the  HF  Benchmark  is  the  development  of  a  generalized 
spatial  smoothing  routine  that  will  accept  arbitrary  ionospheric  profiles.  A  candidate  spatial 
smoothing  algorithm  was  developed  by  Paul  [1991, 1993]  who  demonstrated,  with  satisfying 
results,  a  general  function  that  permitted  fitting  a  set  of  irregularly  spaced  discrete  data.  Fergu¬ 
son  and  Shellman  [1991]  and  Ferguson  [1993]  applied  this  algorithm  to  a  three-dimensional 
array  of  ionospheric  profiles  of  electron  density  with  altitude  distributed  in  a  regular  grid  of 
latitude  and  longitude.  The  fundamental  problem  with  this  technique,  when  applied  to  large  data 
sets,  is  the  large  amount  of  time  required  to  obtain  a  set  of  interpolation  coefficients.  This  limits 
the  flexibility  of  the  resulting  program  because  even  small  changes  in  the  propagation  path 
require  repetition  of  the  whole  iterative  process.  We  now  have  a  new,  more  flexible,  and  faster 
running  technique  employing  a  three-dimensional  cubic  spline  interpolation. 

The  new  procedure  is  quite  straightforward.  The  basic  input  to  the  program  requires 
definition  of  a  reference  path,  usually  between  a  transmitter  and  receiver.  The  program  automati¬ 
cally  sets  up  a  uniform  grid  of  locations  relative  to  this  reference  path.  These  locations  are 
uniformly  spaced  along  the  reference  path  and  along  great  circles  perpendicular  to  it.  The  current 
configuration  is  for  points  along  the  path  to  be  200  km  apart  and  for  the  points  perpendicular  to 
the  path  to  be  at  ±400  and  ±200  km  from  it.  The  electron  densities  are  computed  at  regular 
intervals  of  10  km  in  height.  Following  standard  practice,  it  is  actually  the  logarithm  of  the 
electron  densities  that  is  interpolated. 

Let  us  designate  the  matrix  of  ionospheric  parameters  as/(ij±),  where  i,j,  and  k  represent 
values  along  the  path,  perpendicular  to  the  path,  and  in  height,  respectively.  Since  the  number  of 
points  along  the  path  is  usually  much  greater  than  the  number  of  points  perpendicular  to  the  path 
or  the  number  of  heights,  coefficients  for  the  cubic  spline  interpolation  along  the  reference  path 
(parametric  in  j  and  k)  are  computed  and  stored  only  one  time  for  each  run.  At  an  arbitrary  point 
(xO^z),  we  use  these  coefficients  to  compute  at  set  of  interpolated  values,  representing  the 
original  function  in  the  y-z  plane  at  x.  Then,  by  using  ihestfxijjc)  values,  we  compute  a  set  of 
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coefficients  for  cubic  spline  interpolation  (parametric  in  k)  along  each  of  the  short  paths 
perpendicular  to  the  reference  path.  These  coefficients  are  used  to  interpolate  a  set  of  interpo¬ 
lated  values,  representing  the  function  along  the  z  axis  at  (x,y).  Finally,  we  use  the  fxy{k) 
along  the  z  axis  to  compute  a  set  of  coefficients  for  the  cubic  spline  interpolation  along  the 
z  axis,  which  in  turn,  are  used  to  get  an  interpolated  value,  fxyz  at  the  desired  point.  Since  we  are 
using  cubic  spline  interpolation,  we  also  get  continuous  derivatives  in  all  three  directions 
(necessary  for  the  ray  tracing).  We  achieve  a  lot  of  computational  efficiency  in  solving  for  the 
coefficients  because  the  set  of  equations  which  must  be  solved  are  tri-diagonal,  and  the  coeffi¬ 
cients  of  the  set  of  simultaneous  equations  to  be  solved  for  each  set  of  spline  coefficients  are  the 
same  in  all  three  spatial  directions. 

Four  subroutines  are  used  to  implement  the  interpolation.  The  routine  named  GRD_Model 
takes  the  parameters  of  the  reference  path  and  sets  up  the  matrix  of  ionospheric  data.  This 
routine  also  sets  up  the  interpolation  and  calls  the  interpolation  routine.  The  routine  named 
GRD_Coeff  does  the  cubic  spline  interpolation  for  the  ionospheric  density  and  its  derivatives. 
The  routine  XYThPh  does  coordinate  conversion  between  geographic  coordinates  and  the 
quasi-Cartesian  coordinate  system  used  in  the  interpolation  routine.  This  routine  also  calculates 
the  partial  derivatives  of  geographic  coordinates  with  respect  to  the  quasi-Cartesian  coordinates. 
The  routine  named  PGeoPMag  computes  partial  derivatives  of  geographic  coordinates  with 
respect  to  geomagnetic  coordinates. 

IONOSPHERIC  MODEL 

The  ionospheric  model  used  in  the  program  has  been  updated  to  the  Parameterized  Iono¬ 
spheric  Model  (PIM),  Version  1.05.  (PIM  is  available  from  Computational  Physics  Inc., 

385  Elliot  St.,  Newton,  MA  02164.)  The  realtime  version  of  this  program  is  described  by 
Daniels  [1991].  A  driver  routine,  named  Get_Profile,  has  been  developed  for  this  model.  It  uses 
solar  parameters  such  as  10.7  cm  flux  and  sunspot  number,  supplied  by  the  user,  to  prepare 
inputs  for  the  PIM.  These  user  inputs  are  described  later  in  this  report.  The  driver  routine 
replaces  a  number  of  routines  supplied  as  part  of  the  PIM.  No  changes  other  than  those  required 
for  the  current  operating  system  have  been  made  to  the  subroutines  of  the  PIM  to  allow  quick 
replacement  or  update  of  the  ionospheric  model  in  the  future. 

The  PIM  has  a  few  known  shortcomings.  These  include  a  discontinuity  in  the  E-layer  at  the 
boundary  between  middle  and  high  latitudes.  There  can  be  unusually  steep  density  gradients  and 
bumps  in  the  top  and  bottom  sides  of  the  ionospheric  profile.  These  problems  are  being 
addressed  by  the  developers  of  the  model. 
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RAY  TRACING 


As  in  Ferguson  and  Shellman  [1991]  and  Ferguson  [1993],  the  HF  Benchmark  continues  to 
use  the  Jones  and  Stephenson  [1975]  ray  tracing  model.  This  is  a  versatile  program  with  full 
allowance  for  externally  specified  models  of  the  electron  density,  collision  frequency,  and 
geomagnetic  field.  Every  effort  has  been  made  to  retain  this  flexibility  in  its  implementation  in 
the  HF  Benchmark.  In  particular,  the  input  requirements  are  nearly  identical  to  the  original  Jones 
and  Stephenson  version  with  the  spatial  smoothing  algorithm  being  specified  as  an  ionospheric 
model  named  GRDMODEL.  When  the  program  sees  GRDMODEL  as  the  ionospheric  model, 
the  next  input  record  is  read  as  a  case  identification  consisting  of  one  word,  call  this  word  the 
CASE-ID,  which  forms  the  root  of  the  file  name  for  all  subsequent  outputs.  For  example,  the 
parameters  of  the  ray  paths  are  stored  in  a  file  with  a  name  of  the  form  CASE-ID.RAY  and  the 
execution  log  is  written  to  a  file  names  CASE-ID.LOG.  During  execution  of  the  program,  a 
summary  of  progress  is  written  to  the  screen,  but  all  of  the  detailed  output  is  written  to  this  log 
file.  In  the  present  implementation  of  the  HF  Benchmark  program,  the  calculations  of  the  spatial 
derivatives  of  the  ionospheric  electron  densities  require  that  the  magnetic  field  model  used  in  the 
ray  tracing  be  that  of  a  dipole  (“DIPOLY”  in  the  Jones  and  Stephenson  input  specification). 

In  addition  to  the  original  program  inputs,  a  number  of  additional  parameters,  summarized  in 
figure  1,  have  been  introduced.  The  change  to  the  new  interpolation  algorithm  results  in 
removing  a  number  of  parameters  previously  added  by  Ferguson  and  Shellman  [1991].  The 
value  of  parameter  100  indicates  whether  to  compute  the  ionospheric  profiles  (a  value  of  1)  or  to 
use  an  existing  file  of  ionospheric  profiles  (a  value  of  9).  If  the  value  of  parameter  100  is  1,  then 
the  program  computes  the  required  matrix  of  ionospheric  profiles  and  writes  them  to  a  file 
named  CASE-ID.GRD.  The  meanings  of  parameters  111  through  116  are  self  evident.  Parame¬ 
ter  117  is  a  flag  to  indicate  the  direction  of  the  solar  magnetic  field.  By.  Parameters  118  and  119 
set  lower  and  upper  limits,  respectively,  on  the  heights  of  the  ionosphere.  If  the  value  of 
parameter  100  is  9,  the  values  of  the  parameters  118  and  119  may  be  different  from  those  used  to 
compute  the  ionospheric  profiles  found  in  the  file.  This  allows  the  user  to  experiment  with 
implicit  removal  of  E  layer  effects  or  to  lower  the  effective  top  of  the  ionosphere. 

Parameters  130, 131  and  132  set  up  the  parameters  of  the  ionospheric  display  for  the 
reference  path.  The  first  parameter  130  simply  initiates  the  block  of  input.  Parameter  131  defines 
the  quantity  to  be  plotted.  Two  views  of  this  quantity  are  produced:  the  first  view  is  in  a  vertical 
plane  passing  through  the  endpoints  of  the  reference  path;  the  second  is  in  the  horizontal  plane  at 
a  height  defined  by  parameter  132.  The  original  intention  was  to  have  the  HF  Benchmark 
program  generate  this  display,  but  it  is  much  easier  to  have  a  separate  program,  named  GridPlot, 
do  this  for  now. 

Because  the  ionosphere  is  defined  over  a  limited  area,  the  ray  tracing  routines  have  been 
modified  to  test  for  the  ray  path  moving  outside  of  the  defined  area.  An  appropriate  message  is 
printed  into  the  log  file  when  this  occurs.  Occasionally,  if  the  index  of  refraction  becomes 
negative  and  the  problem  carmot  be  circumvented,  we  terminate  the  calculations  for  the  current 
ray  and  print  an  appropriate  message  into  the  log  file.  A  number  of  other  minor  procedural 
modifications  have  been  made  to  the  ray  tracing  program  to  facilitate  miming  it  on  a  PC. 
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Figure  1.  Parameters  added  to  the  ray  trace  input  array  to  support  the  HF  Benchmark. 


SAMPLE  PROBLEM 

A  sample  program  is  presented  here  to  illustrate  the  unification  of  the  individual  models.  The 
input  parameters  for  the  program  are  shown  in  figure  2.  The  path  is  defined  in  geographic 
coordinates  from  68.5°N,  32.5°E  to  55.5®N,  117.1  °W  (parameters  4, 5, 19,  and  20).  The 
computed  geogri^hic  bearing  angle  from  the  transmitter  to  the  receiver  is  339.3°  and  the  path 
length  is  6000  km;  therefore,  the  initial  azimuth  of  the  rays  is  set  to  340°  clockvose  from  North 
(parameter  11).  Ordinary  rays  will  be  laimched  at  this  azimuth  at  elevation  angles  from  0°  to 
15°  (parameters  15  and  16)  in  steps  of  3°  (parameter  17).  The  ionospheric  profiles  are  generated 
for  15  December  1991  at  1500  UT  (parameters  11  through  114).  The  10  ran  flux  is  70,  the  Up 
is  2,  and  the  direction  of  the  sun’s  magnetic  field.  By,  is  positive  (parameters  115  through  117). 
Electron  densities  are  computed  from  90  to  400  km  (parameters  118  and  119).  If  so  requested, 
plots  of  the  ionospheric  parameters  will  give  plasma  frequency  in  the  vertical  plane  and  in  a 
horizontal  plane  at  200  km  (parameters  131  and  132).  Plots  of  the  iononspheric  parameter  and  of 
t:he  ray  paths  will  be  made  in  the  plane  of  the  reference  path  (parameters  82  through  85)  and  the 
horizontal  displays  will  have  an  extent  of  300  km  (parameter  87). 

Contour  plots  of  the  plasma  frequency  in  the  vertical  plane  along  the  path  are  shown  in  the 
top  panel  of  figure  3,  and  contours  of  die  plasma  frequency  in  the  horizontal  plane  at  an  altitude 
of  200  km  are  shown  in  the  bottom  panel.  This  display  is  generated  by  an  auxiliary  program 
named  GridPlot  This  program  uses  the  same  input  as  the  HF  Benchmark  and  calculates  the 
selected  ionospheric  parametra  by  using  the  file  containing  the  ionospheric  profiles.  It  is  clear 
that  the  ionosphere  is  changing  in  all  three  dimensions.  A  new  feature  of  the  gri^hical  ouqiut  is 
the  use  of  filled  contour  plots  instead  of  the  contour  lines  used  in  the  past.  In  this  gray  shade 
version,  there  is  not  enough  resolution  to  visualize  the  complexity  of  the  ionospheric  variation; 
however,  this  complexity  is  evident  when  the  displays  are  viewed  in  color. 
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Figure  2.  Inputs  for  sample  problem. 
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Figure  3.  Plots  of  contours  of  constant  values  of  the  logarithm  of  the  plasma  frequency  in  the 
vertical  plane  along  the  great  circle  path  passing  through  the  transmitter  and  receive*  in  the  hori¬ 
zontal  plane  at  200  km. 
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Ray  paths  fcM*  the  extraordinary  ray  at  4  MHz  are  shown  in  figure  4.  The  top  panel  shows  the 
projection  of  the  ray  paths  onto  the  vertical  plane  passing  through  the  transmitter  and  the 
receiver.  The  bottom  panel  of  figure  4  shows  the  projection  of  the  ray  paths  onto  the  ground  with 
locations  above  the  horizontal  axis  being  lo  the  west  of  the  propagation  path.  We  see  that  the 
off-path  variation  of  the  ionosphere  causes  some  deflection  of  the  ray  paths. 


Figure  4.  Plots  of  the  ray  paths  for  the  extraordinary  ray  at  4  MHz  in  the  vertical  plane  passing 
through  the  transmitter  and  the  receiver  (top  panel)  and  in  the  horizontal  plane  (bottom  panel). 
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One  of  the  displays  is  a  global  view  showing  contours  of  the  ionospheric  profile  in  a  fixed 
meridian,  and  the  other  is  a  view  along  short  paths  (a  subset  of  the  global  view).  In  each  display, 
contours  of  the  electron  density  taken  directly  from  the  ionospheric  model  are  shown.  An 
example  of  the  global  display  is  shown  in  figure  6.  The  map  projection  is  orthographic  with  the 
observation  point  in  the  center  of  the  display.  The  location  of  this  point  is  specified  by  the  user. 
Superimposed  on  the  circle  representing  the  earth  are  a  series  of  lines  indicating  land  masses. 
There  are  a  series  of  lines,  ranging  fi-om  solid,  through  dashes  of  various  lengths,  to  dotted. 

These  lines  indicate  the  position  of  solar  zenith  angles,  ranging  from  90  degrees  (solid)  to 
99  degrees  (dotted),  to  show  where  the  day-night  transition  is  located.  If  it  is  visible  from  the 
observation  side  of  the  globe,  the  sub-solar  point  is  also  shown.  The  vertical  scale  for  points 
above  the  earth’s  surface  is  specified  by  the  user  setting  the  number  of  inches  between  the  earth’s 
surface  and  the  top  of  the  model  ionosphere.  The  auroral  bulges  are  visible  and  the  differences 
between  the  daytime  and  nighttime  ionospheres  are  evident.  The  data  set  required  to  make  this 
display  is  generated  by  a  program  named  Globe,  and  the  plots  are  generated  by  a  program  named 
GlobePlot. 

The  second  display  is  a  modification  of  the  first,  a  sample  of  which  is  shown  in  figure  7.  The 
primary  difference  from  the  global  display  is  the  shorter  path  length  which  allows  scaling  of  the 
vertical  dimension  closer  to  the  horizontal  dimension.  Consistent  with  the  names  used  above,  the 
data  set  required  to  make  this  display  is  generated  by  a  program  named  Cntr,  and  the  plots  are 
generated  by  a  program  named  CntrPlot.  This  display  is  similar  to  the  one  generated  by  the 
HF  Benchmark  program,  lacking  the  horizontal  display;  furthermore,  the  Benchmark  ^play  is 
the  result  of  the  spline  interpolation  used  in  the  ray  tracing,  while  the  display  shown  in  figure  7 
is  the  result  of  computing  ionospheric  profiles  closely  spaced  in  distance. 

FUTURE  PLANS 


Our  primary  goal  in  the  next  year  is  to  exercise  the  program  over  propagation  paths  we  have 
measured  so  that  we  can  make  comparisons.  We  can  then  begin  to  sort  out  what  is  important  in 
the  ionospheric  model  and  what  is  not  relative  to  HP  propagation.  There  are  known  shortcom¬ 
ings  in  the  ionospheric  model,  but  we  must  wait  for  its  developers  to  solve  them.  A  minor 
problem  in  the  ray  tracing  must  be  fixed,  namely,  the  ray  path  is  often  extrapolated  during  the 
calculation.  Sometimes  this  extr{q)olation  puts  a  calculation  point  outside  of  the  defined 
ionosphere  and  we  get  a  negative  index  of  refraction.  Currently,  we  simply  stop  tracing  that 
particular  ray  and  continue  with  the  next  one.  We  have  to  expand  the  checks  on  the  ray  being  out 
of  bounds  to  solve  this  problem.  The  current  problem  is  quite  simple  to  use  and  the  computa¬ 
tional  throughput  is  acceptable  for  small  scale  studies  of  HF  propagation. 
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